Load Data

dataset <- read.delim("raw_data/FigureS5L.txt", stringsAsFactors = FALSE)

dataset$genotype <- gsub(" ", "", dataset$genotype )
dataset$genotype <- factor(dataset$genotype)
dataset$Experiment <- factor(rep(paste0("exp", 1:(length(dataset$genotype)/length(levels(dataset$genotype)))), each=length(unique(dataset$genotype))))

dataset$UID <- factor(paste(dataset$Experiment, dataset$genotype))

# wide format
kable(dataset, row.names = F)
genotype NT bleomycin_5uM bleomycin_10uM bleomycin_50uM Experiment UID
WT 1160 810 567 151 exp1 exp1 WT
YFP-ALC1 1005 792 585 192 exp1 exp1 YFP-ALC1
WT 1056 838 378 258 exp2 exp2 WT
YFP-ALC1 2132 1742 983 399 exp2 exp2 YFP-ALC1
WT 1604 1036 558 184 exp3 exp3 WT
YFP-ALC1 1714 1325 985 311 exp3 exp3 YFP-ALC1
library(reshape2)
# reshape to long format
dataset <- melt(dataset, variable.name = "Treatment", value.name = "Counts")

dataset$genotype <- relevel(dataset$genotype, ref = "WT")
dataset$UID <- relevel(dataset$UID, ref = "exp1 WT")

dataset$Bleomycin <- gsub("NT","1",dataset$Treatment)
dataset$Bleomycin <- gsub("bleomycin_|uM","",dataset$Bleomycin)
dataset$Bleomycin <- log10(as.integer(dataset$Bleomycin))




dataset$Offset <- NA
for(uid in levels(dataset$UID)){
        dataset$Offset[dataset$UID == uid] <- mean(dataset$Counts[dataset$UID == uid])
}

dataset$NormCounts <- dataset$Counts / dataset$Offset



dataset$Offset2 <- NA
for(gid in levels(dataset$genotype)){
        dataset$Offset2[dataset$genotype == gid] <- mean(dataset$NormCounts[dataset$genotype == gid & dataset$Bleomycin == 0])
}

dataset$NormCounts2 <- dataset$NormCounts / dataset$Offset2



# long format
kable(dataset, row.names = F)
genotype Experiment UID Treatment Counts Bleomycin Offset NormCounts Offset2 NormCounts2
WT exp1 exp1 WT NT 1160 0.00000 672.00 1.7261905 1.764286 0.9784074
YFP-ALC1 exp1 exp1 YFP-ALC1 NT 1005 0.00000 643.50 1.5617716 1.588615 0.9831029
WT exp2 exp2 WT NT 1056 0.00000 632.50 1.6695652 1.764286 0.9463121
YFP-ALC1 exp2 exp2 YFP-ALC1 NT 2132 0.00000 1314.00 1.6225266 1.588615 1.0213469
WT exp3 exp3 WT NT 1604 0.00000 845.50 1.8971023 1.764286 1.0752805
YFP-ALC1 exp3 exp3 YFP-ALC1 NT 1714 0.00000 1083.75 1.5815456 1.588615 0.9955502
WT exp1 exp1 WT bleomycin_5uM 810 0.69897 672.00 1.2053571 1.764286 0.6831983
YFP-ALC1 exp1 exp1 YFP-ALC1 bleomycin_5uM 792 0.69897 643.50 1.2307692 1.588615 0.7747438
WT exp2 exp2 WT bleomycin_5uM 838 0.69897 632.50 1.3249012 1.764286 0.7509560
YFP-ALC1 exp2 exp2 YFP-ALC1 bleomycin_5uM 1742 0.69897 1314.00 1.3257230 1.588615 0.8345152
WT exp3 exp3 WT bleomycin_5uM 1036 0.69897 845.50 1.2253105 1.764286 0.6945078
YFP-ALC1 exp3 exp3 YFP-ALC1 bleomycin_5uM 1325 0.69897 1083.75 1.2226067 1.588615 0.7696056
WT exp1 exp1 WT bleomycin_10uM 567 1.00000 672.00 0.8437500 1.764286 0.4782388
YFP-ALC1 exp1 exp1 YFP-ALC1 bleomycin_10uM 585 1.00000 643.50 0.9090909 1.588615 0.5722539
WT exp2 exp2 WT bleomycin_10uM 378 1.00000 632.50 0.5976285 1.764286 0.3387367
YFP-ALC1 exp2 exp2 YFP-ALC1 bleomycin_10uM 983 1.00000 1314.00 0.7480974 1.588615 0.4709118
WT exp3 exp3 WT bleomycin_10uM 558 1.00000 845.50 0.6599645 1.764286 0.3740689
YFP-ALC1 exp3 exp3 YFP-ALC1 bleomycin_10uM 985 1.00000 1083.75 0.9088812 1.588615 0.5721219
WT exp1 exp1 WT bleomycin_50uM 151 1.69897 672.00 0.2247024 1.764286 0.1273617
YFP-ALC1 exp1 exp1 YFP-ALC1 bleomycin_50uM 192 1.69897 643.50 0.2983683 1.588615 0.1878167
WT exp2 exp2 WT bleomycin_50uM 258 1.69897 632.50 0.4079051 1.764286 0.2312013
YFP-ALC1 exp2 exp2 YFP-ALC1 bleomycin_50uM 399 1.69897 1314.00 0.3036530 1.588615 0.1911433
WT exp3 exp3 WT bleomycin_50uM 184 1.69897 845.50 0.2176227 1.764286 0.1233489
YFP-ALC1 exp3 exp3 YFP-ALC1 bleomycin_50uM 311 1.69897 1083.75 0.2869666 1.588615 0.1806395

Plot Data

library(ggplot2)

# raw data
ggplot(dataset, aes(x=Bleomycin, y=Counts)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE, aes(colour=genotype)) +
        geom_point(aes(colour=genotype, shape=Experiment), size=2) +        
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_shape_manual(values=15:20) +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)")+
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)")+
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

library(Cairo)

cairo_pdf("FigureS5L.pdf", width = 5, height = 4, family = "Arial")

ggplot(dataset, aes(x=Bleomycin, y=NormCounts2)) + 
        theme_bw() +
        theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(), 
              axis.line = element_line(colour = "black"), text = element_text(size=14),
              panel.border = element_blank(), panel.background = element_blank()) +
        geom_point(aes(colour = genotype)) +
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=TRUE, aes(colour = genotype), fill='#CCCCCC') +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        ylab(label = "Normalized Counts") +
        scale_color_manual(values=c("#000000","#808000"))

dev.off()
## quartz_off_screen 
##                 2

Models

library(MASS)
library(DHARMa)
library(lme4)
library(lmerTest)
library(bbmle)

Linear formula

fit1 <- lm(Counts ~ Experiment + Bleomycin*genotype, data = dataset)
print(summary(fit1))
## 
## Call:
## lm(formula = Counts ~ Experiment + Bleomycin * genotype, data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -476.82 -107.86  -15.27  130.19  500.52 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1064.5      141.8   7.509 5.97e-07 ***
## Experimentexp2                315.5      123.7   2.550   0.0201 *  
## Experimentexp3                306.9      123.7   2.481   0.0232 *  
## Bleomycin                    -653.7      117.1  -5.583 2.68e-05 ***
## genotypeYFP-ALC1              417.3      173.2   2.410   0.0269 *  
## Bleomycin:genotypeYFP-ALC1   -141.6      165.6  -0.855   0.4038    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 247.4 on 18 degrees of freedom
## Multiple R-squared:  0.8398, Adjusted R-squared:  0.7953 
## F-statistic: 18.88 on 5 and 18 DF,  p-value: 1.355e-06
cat("AIC: ", AIC(fit1))
## AIC:  339.7383
simres <- simulateResiduals(fittedModel = fit1)
plot(simres)

fit2 <- lm(NormCounts ~ Bleomycin*genotype, data = dataset)
print(summary(fit2))
## 
## Call:
## lm(formula = NormCounts ~ Bleomycin * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.26678 -0.05879 -0.01889  0.09365  0.20857 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 1.76526    0.06366  27.729  < 2e-16 ***
## Bleomycin                  -0.90086    0.06087 -14.799 3.08e-12 ***
## genotypeYFP-ALC1           -0.10406    0.09003  -1.156    0.261    
## Bleomycin:genotypeYFP-ALC1  0.12249    0.08608   1.423    0.170    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1286 on 20 degrees of freedom
## Multiple R-squared:  0.9503, Adjusted R-squared:  0.9429 
## F-statistic: 127.5 on 3 and 20 DF,  p-value: 3.315e-13
cat("AIC: ", AIC(fit2))
## AIC:  -24.70436
simres <- simulateResiduals(fittedModel = fit2)
plot(simres)

fit3 <- lm(NormCounts2 ~ Bleomycin*genotype, data = dataset)
print(summary(fit3))
## 
## Call:
## lm(formula = NormCounts2 ~ Bleomycin * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.15121 -0.03700 -0.01070  0.05473  0.13129 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 1.00055    0.03746  26.709  < 2e-16 ***
## Bleomycin                  -0.51061    0.03582 -14.255 6.13e-12 ***
## genotypeYFP-ALC1            0.04514    0.05298   0.852    0.404    
## Bleomycin:genotypeYFP-ALC1  0.02064    0.05066   0.408    0.688    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07569 on 20 degrees of freedom
## Multiple R-squared:  0.9517, Adjusted R-squared:  0.9445 
## F-statistic: 131.5 on 3 and 20 DF,  p-value: 2.479e-13
cat("AIC: ", AIC(fit3))
## AIC:  -50.15771
simres <- simulateResiduals(fittedModel = fit3)
plot(simres)

fit4 <- lmer(Counts ~ Bleomycin*genotype + (1|UID), data = dataset)
print(summary(fit4))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ Bleomycin * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 282.7
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.91400 -0.53151 -0.01178  0.60395  1.83213 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 55274    235.1   
##  Residual             36717    191.6   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                            Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)                1271.943    165.583    6.393   7.682 0.000187 ***
## Bleomycin                  -653.663     90.676   16.000  -7.209 2.08e-06 ***
## genotypeYFP-ALC1            417.333    234.169    6.393   1.782 0.121961    
## Bleomycin:genotypeYFP-ALC1 -141.556    128.235   16.000  -1.104 0.285973    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Blmycn gYFP-A
## Bleomycin   -0.465              
## gntYFP-ALC1 -0.707  0.329       
## Bl:YFP-ALC1  0.329 -0.707 -0.465
cat("AIC: ", AIC(fit4))
## AIC:  294.7363
simres <- simulateResiduals(fittedModel = fit4)
plot(simres)

Quadratic formula

fit5 <- lm(Counts ~ Experiment + poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit5))
## 
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 2) * genotype, 
##     data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -421.90 -106.95  -25.02  149.03  445.60 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            509.21     105.15   4.842  0.00018 ***
## Experimentexp2                         315.50     128.79   2.450  0.02618 *  
## Experimentexp3                         306.88     128.79   2.383  0.02993 *  
## poly(Bleomycin, 2)1                  -1953.50     364.27  -5.363 6.34e-05 ***
## poly(Bleomycin, 2)2                     92.26     364.27   0.253  0.80327    
## genotypeYFP-ALC1                       297.08     105.15   2.825  0.01219 *  
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  -423.04     515.15  -0.821  0.42360    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  -361.30     515.15  -0.701  0.49316    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 257.6 on 16 degrees of freedom
## Multiple R-squared:  0.8457, Adjusted R-squared:  0.7782 
## F-statistic: 12.53 on 7 and 16 DF,  p-value: 1.973e-05
cat("AIC: ", AIC(fit5))
## AIC:  342.8408
simres <- simulateResiduals(fittedModel = fit5)
plot(simres)

fit6 <- lm(NormCounts ~ poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit6))
## 
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 2) * genotype, data = dataset)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.242931 -0.041198  0.004446  0.065468  0.213157 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           1.000e+00  3.629e-02  27.559 3.58e-16 ***
## poly(Bleomycin, 2)1                  -2.692e+00  1.778e-01 -15.145 1.10e-11 ***
## poly(Bleomycin, 2)2                   1.168e-01  1.778e-01   0.657    0.519    
## genotypeYFP-ALC1                      4.392e-17  5.132e-02   0.000    1.000    
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  3.661e-01  2.514e-01   1.456    0.163    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1 -3.987e-01  2.514e-01  -1.586    0.130    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1257 on 18 degrees of freedom
## Multiple R-squared:  0.9573, Adjusted R-squared:  0.9454 
## F-statistic: 80.71 on 5 and 18 DF,  p-value: 1.097e-11
cat("AIC: ", AIC(fit6))
## AIC:  -24.34165
simres <- simulateResiduals(fittedModel = fit6)
plot(simres)

fit7 <- lm(NormCounts2 ~ poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit7))
## 
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 2) * genotype, data = dataset)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.137694 -0.024008  0.002698  0.039744  0.120818 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           0.56680    0.02115  26.794 5.88e-16 ***
## poly(Bleomycin, 2)1                  -1.52597    0.10363 -14.724 1.76e-11 ***
## poly(Bleomycin, 2)2                   0.06622    0.10363   0.639   0.5309    
## genotypeYFP-ALC1                      0.06268    0.02992   2.095   0.0506 .  
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  0.06170    0.14656   0.421   0.6788    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1 -0.24363    0.14656  -1.662   0.1138    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07328 on 18 degrees of freedom
## Multiple R-squared:  0.9593, Adjusted R-squared:  0.948 
## F-statistic: 84.83 on 5 and 18 DF,  p-value: 7.161e-12
cat("AIC: ", AIC(fit7))
## AIC:  -50.24129
simres <- simulateResiduals(fittedModel = fit7)
plot(simres)

fit8 <- lmer(Counts ~ poly(Bleomycin,2)*genotype + (1|UID), data = dataset)
print(summary(fit8))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 2) * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 251.1
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.59469 -0.57627 -0.04532  0.59448  1.51209 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 54685    233.8   
##  Residual             39074    197.7   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                                      Estimate Std. Error       df t value
## (Intercept)                            716.67     146.57     4.00   4.889
## poly(Bleomycin, 2)1                  -1953.49     279.55    14.00  -6.988
## poly(Bleomycin, 2)2                     92.26     279.55    14.00   0.330
## genotypeYFP-ALC1                       297.08     207.29     4.00   1.433
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  -423.05     395.34    14.00  -1.070
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  -361.30     395.34    14.00  -0.914
##                                      Pr(>|t|)    
## (Intercept)                           0.00811 ** 
## poly(Bleomycin, 2)1                  6.37e-06 ***
## poly(Bleomycin, 2)2                   0.74625    
## genotypeYFP-ALC1                      0.22509    
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  0.30269    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  0.37624    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) pl(B,2)1 pl(B,2)2 gYFP-A p(B,2)1:
## ply(Blm,2)1  0.000                                  
## ply(Blm,2)2  0.000  0.000                           
## gntYFP-ALC1 -0.707  0.000    0.000                  
## p(B,2)1:YFP  0.000 -0.707    0.000    0.000         
## p(B,2)2:YFP  0.000  0.000   -0.707    0.000  0.000
cat("AIC: ", AIC(fit8))
## AIC:  267.1439
simres <- simulateResiduals(fittedModel = fit8)
plot(simres)

Cubic formula

fit9 <- lm(Counts ~ Experiment + poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit9))
## 
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 3) * genotype, 
##     data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -404.54  -95.08   -6.06  104.79  406.96 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            509.208    105.905   4.808 0.000278 ***
## Experimentexp2                         315.500    129.707   2.432 0.029009 *  
## Experimentexp3                         306.875    129.707   2.366 0.032948 *  
## poly(Bleomycin, 3)1                  -1953.495    366.867  -5.325 0.000107 ***
## poly(Bleomycin, 3)2                     92.264    366.867   0.251 0.805088    
## poly(Bleomycin, 3)3                    346.343    366.867   0.944 0.361149    
## genotypeYFP-ALC1                       297.083    105.905   2.805 0.014036 *  
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  -423.045    518.828  -0.815 0.428506    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  -361.300    518.828  -0.696 0.497592    
## poly(Bleomycin, 3)3:genotypeYFP-ALC1    -1.664    518.828  -0.003 0.997487    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 259.4 on 14 degrees of freedom
## Multiple R-squared:  0.8631, Adjusted R-squared:  0.775 
## F-statistic: 9.804 on 9 and 14 DF,  p-value: 0.0001174
cat("AIC: ", AIC(fit9))
## AIC:  343.9775
simres <- simulateResiduals(fittedModel = fit9)
plot(simres)

fit10 <- lm(NormCounts ~ poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit10))
## 
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 3) * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.10726 -0.04199 -0.01795  0.05358  0.14330 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           1.000e+00  2.481e-02  40.307  < 2e-16 ***
## poly(Bleomycin, 3)1                  -2.692e+00  1.215e-01 -22.151 1.97e-13 ***
## poly(Bleomycin, 3)2                   1.168e-01  1.215e-01   0.961 0.350736    
## poly(Bleomycin, 3)3                   4.929e-01  1.215e-01   4.056 0.000918 ***
## genotypeYFP-ALC1                      8.710e-17  3.509e-02   0.000 1.000000    
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  3.661e-01  1.719e-01   2.130 0.049063 *  
## poly(Bleomycin, 3)2:genotypeYFP-ALC1 -3.987e-01  1.719e-01  -2.319 0.033931 *  
## poly(Bleomycin, 3)3:genotypeYFP-ALC1 -1.938e-01  1.719e-01  -1.128 0.276093    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.08594 on 16 degrees of freedom
## Multiple R-squared:  0.9823, Adjusted R-squared:  0.9745 
## F-statistic: 126.5 on 7 and 16 DF,  p-value: 8.384e-13
cat("AIC: ", AIC(fit10))
## AIC:  -41.41727
simres <- simulateResiduals(fittedModel = fit10)
plot(simres)

fit11 <- lm(NormCounts2 ~ poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit11))
## 
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 3) * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.06752 -0.02410 -0.01047  0.03373  0.08122 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           0.56680    0.01442  39.314  < 2e-16 ***
## poly(Bleomycin, 3)1                  -1.52597    0.07063 -21.605  2.9e-13 ***
## poly(Bleomycin, 3)2                   0.06622    0.07063   0.938  0.36241    
## poly(Bleomycin, 3)3                   0.27939    0.07063   3.956  0.00113 ** 
## genotypeYFP-ALC1                      0.06268    0.02039   3.074  0.00726 ** 
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  0.06170    0.09989   0.618  0.54549    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1 -0.24363    0.09989  -2.439  0.02675 *  
## poly(Bleomycin, 3)3:genotypeYFP-ALC1 -0.09112    0.09989  -0.912  0.37520    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04994 on 16 degrees of freedom
## Multiple R-squared:  0.9832, Adjusted R-squared:  0.9758 
## F-statistic: 133.7 on 7 and 16 DF,  p-value: 5.446e-13
cat("AIC: ", AIC(fit11))
## AIC:  -67.4721
simres <- simulateResiduals(fittedModel = fit11)
plot(simres)

fit12 <- lmer(Counts ~ poly(Bleomycin,3)*genotype + (1|UID), data = dataset)
print(summary(fit12))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 3) * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 221.8
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.55170 -0.46246 -0.07901  0.43109  1.35742 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 55544    235.7   
##  Residual             35638    188.8   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                                       Estimate Std. Error        df t value
## (Intercept)                            716.667    146.575     4.000   4.889
## poly(Bleomycin, 3)1                  -1953.495    266.975    12.000  -7.317
## poly(Bleomycin, 3)2                     92.264    266.975    12.000   0.346
## poly(Bleomycin, 3)3                    346.343    266.975    12.000   1.297
## genotypeYFP-ALC1                       297.083    207.289     4.000   1.433
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  -423.045    377.560    12.000  -1.120
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  -361.300    377.560    12.000  -0.957
## poly(Bleomycin, 3)3:genotypeYFP-ALC1    -1.664    377.560    12.000  -0.004
##                                      Pr(>|t|)    
## (Intercept)                           0.00811 ** 
## poly(Bleomycin, 3)1                  9.26e-06 ***
## poly(Bleomycin, 3)2                   0.73563    
## poly(Bleomycin, 3)3                   0.21892    
## genotypeYFP-ALC1                      0.22509    
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  0.28444    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  0.35748    
## poly(Bleomycin, 3)3:genotypeYFP-ALC1  0.99656    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) pl(B,3)1 pl(B,3)2 pl(B,3)3 gYFP-A p(B,3)1: p(B,3)2:
## ply(Blm,3)1  0.000                                                    
## ply(Blm,3)2  0.000  0.000                                             
## ply(Blm,3)3  0.000  0.000    0.000                                    
## gntYFP-ALC1 -0.707  0.000    0.000    0.000                           
## p(B,3)1:YFP  0.000 -0.707    0.000    0.000    0.000                  
## p(B,3)2:YFP  0.000  0.000   -0.707    0.000    0.000  0.000           
## p(B,3)3:YFP  0.000  0.000    0.000   -0.707    0.000  0.000    0.000
cat("AIC: ", AIC(fit12))
## AIC:  241.8309
simres <- simulateResiduals(fittedModel = fit12)
plot(simres)

Compare Results

ICtab(fit1,fit2,fit3,fit4,
      fit5,fit6,fit7,fit8,
      fit9,fit10,fit11,fit12,
      base=T)
##       AIC   dAIC  df
## fit11 -67.5   0.0 9 
## fit7  -50.2  17.2 7 
## fit3  -50.2  17.3 5 
## fit10 -41.4  26.1 9 
## fit2  -24.7  42.8 5 
## fit6  -24.3  43.1 7 
## fit12 241.8 309.3 10
## fit8  267.1 334.6 8 
## fit4  294.7 362.2 6 
## fit1  339.7 407.2 7 
## fit5  342.8 410.3 9 
## fit9  344.0 411.4 11

Final Result

fit <- fit11

output <- coef(summary(fit))
output <- output[grep("Bleomycin", rownames(output)),]

rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype",  paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1],  sep = " in " )


# suggested result table
kable(output, row.names = T)
Estimate Std. Error t value Pr(>|t|)
Bleomycin1 in WT -1.5259664 0.0706299 -21.6050932 0.0000000
Bleomycin2 in WT 0.0662200 0.0706299 0.9375630 0.3624070
Bleomycin3 in WT 0.2793890 0.0706299 3.9556733 0.0011329
Bleomycin1: WT vs. YFP-ALC1 0.0616959 0.0998858 0.6176642 0.5454882
Bleomycin2: WT vs. YFP-ALC1 -0.2436254 0.0998858 -2.4390393 0.0267529
Bleomycin3: WT vs. YFP-ALC1 -0.0911161 0.0998858 -0.9122025 0.3752023
write.table(output, file = "FigureS5L_Stats.txt", quote = F, sep = "\t", row.names = T, col.names = NA)

Anova

fit11a <- lm(NormCounts2 ~ poly(Bleomycin,3)*genotype, data = dataset)
fit11b <- lm(NormCounts2 ~ poly(Bleomycin,3)+genotype, data = dataset)

# anova table
anova(fit11a, fit11b)
## Analysis of Variance Table
## 
## Model 1: NormCounts2 ~ poly(Bleomycin, 3) * genotype
## Model 2: NormCounts2 ~ poly(Bleomycin, 3) + genotype
##   Res.Df      RSS Df Sum of Sq      F Pr(>F)
## 1     16 0.039909                           
## 2     19 0.057774 -3 -0.017865 2.3875 0.1071